Emergence of self-sustained patterns in small-world excitable media 
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Motivated by recent observations that long-range connections (LRC) play a role in various brain 
phenomena, we have observed two distinct dynamical transitions in the activity of excitable media 
where waves propagate both between neighboring regions and through LRC. When LRC density p 
is low, single or multiple spiral waves are seen to emerge and cover the entire system. This state is 
self-sustaining and robust against perturbations. At p = p l c the spirals are suppressed and there is 
a transition to a spatially homogeneous, temporally periodic state. Finally, above p = p™, activity 
ceases after a brief transient. 

PACS numbers: 05.65.+b, 87.18.Bb, 05.45.Jn, 89.75.Hc 



Pattern formation in excitable media has recently been 
a very exciting area of research, not least because it is 
observed in a wide variety of natural systems, ranging 
from spiral waves in mammalian brain [l[ to chemical 
systems such as the Belousov-Zhabotinsky reaction [2j. 
The functional role of such patterns in biological systems 
makes it imperative to understand better the conditions 
in which they can spontaneously emerge. For example, 
these patterns have been implicated in the genesis of life- 
threatening arrhythmias in the heart [3(], while in the 
brain they have been thought to provide a spatial frame- 
work for cortical oscillations Until now, work has 
mostly focussed on generating patterns in excitable me- 
dia by using stochastic stimulation [4]. In this case, how- 
ever, the activity is not really self-sustained as the noise is 
akin to external intervention necessary for the initiation 
and persistence of spiral waves [1, EJ. If such patterns 
are to be seen as spontaneously emerging from arbitrary 
initial conditions, then the pattern formation should be 
an outcome of the internal structure of the system alone. 
Furthermore, small variations in this structure may re- 
sult in transitions between different dynamical regimes 
characterized by distinct spatiotemporal patterns. 



Here we have considered excitable systems which have 
a regular topology with cells communicating only with 
nearest neighbors, but where there are a few random 
long-range connections, linking spatially distant cells. 
Such "small-world" topologies have been observed in a 
large number of real world systems [7|], and have also 
been associated with self-sustained activity in a chain of 
model neurons [8[ as well as periodic epidemic patterns in 
disease spreading [9]. However, while it has been shown 
that such long-range connections do play a role in main- 
taining spiral waves that already exist [10], so far there 
have been very few attempts at showing the emergence of, 
and transitions between, different types of self-sustained 
patterns as a result of such topology. The observation of 
spontaneous pattern formation in natural systems where 
sparse long-range connections coexist with fairly regu- 
lar underlying connection topology, like in the brain [ll| , 



points towards intriguing possibilities for the functional 
role of these connections. In this Letter we have investi- 
gated a generic model of excitable media with increasing 
density of random long-range connections, and we report 
the existence of two qualitatively different regimes of self- 
sustained pattern formation. The correspondence of the 
observed patterns with those observed in nature, e.g., 
epileptic bursts and seizures [12[, as well as their depen- 
dence on the topological structure of connections, makes 
these results highly relevant in our view. 

The model we consider here consists of a two- 
dimensional array of N x TV excitable cells coupled diffu- 
sively, 
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Here D is the diffusion coefficient and the dynamics of 
individual cells are described by a pair of variables x t: yt, 
evolving according to a discrete-time model of generic 
excitable media 11311: 
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where we have fixed the parameters as a = 0.89, b = 
0.6, c = 0.28, and k — 0.02. When this system is ex- 
cited with a suprathreshold stimulation, the fast variable 
x shows an abrupt increase. This triggers changes in the 
slow variable y, such that the state of the cell is grad- 
ually brought down to that of the resting state. Once 
excited, the cell remains impervious to stimulation up to 
a refractory period, the duration of which is governed by 
the parameter a. Neighboring cells communicate exci- 
tation to each other with a strength proportional to D, 
chosen D = 0.2 for most of our simulations. In addition 
to the diffusive coupling, we introduce long-range con- 
nections such that each cell receives a connection from 
a randomly chosen cell with probability p. The strength 
of this connection is chosen to be the same as that of 
the nearest neighbors, i.e., D/4. These random long- 
range links can be either quenched (i.e., chosen initially 
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FIG. 1: Emergence of spiral waves, recorded in a simulation 
of a N 2 = 128 x 128 system, with shortcut probability p — 
0.25. After initialization, the dynamics of the system can be 
characterized by circular waves. At around t ~ 1,500 time 
steps, a spiral wave is spontaneously created and subsequently 
takes over the dynamics. Panels (a) and (b) show the state of 
the system at different times in terms of the fast variable x. 
Colors indicate excitation level of the cells, (c) Time series of 
the average activity, i.e. the fraction of cells with x > 0.9. 

and kept fixed for the duration of the simulation), or an- 
nealed (i.e., randomly created at each time step). While 
the results reported below are for annealed random links, 
we could observe no qualitative difference between these 
two cases. Note that we have used both periodic and 
absorbing boundary conditions for the system, and ob- 
served no significant difference in the results. All the 
following results were obtained for absorbing boundary 
conditions. 

In the simulations, the observed patterns formed in 
the system for varying values of the shortcut density p 
can, generally speaking, be divided into three categories. 
First, below the lower critical probability, i.e., < p < p l c: 
the state of the system after an initial transient period is 
characterized by self-sustaining single or multiple spiral 
waves covering the entire system. Second, at p = p l c , the 
spiral wave mode is suppressed and the system under- 
goes a transition to a regime in which a large fraction of 
the system gets simultaneously active, and then refrac- 
tory, in a periodic manner. Third, when the value of p 
is increased above a system-size-dependent upper critical 
probability^, the self-sustained activity ceases and the 
system falls into the absorbing state where x 1 ^ = 0, Vi, j. 

In the beginning of each simulation run the system was 
initialized such that x — 1 for a small number of cells, 
and x = for the rest. For small shortcut probabilities 
p <C p l c , upon starting the simulation, multiple coexist- 
ing circular excitation waves were seen to emerge (see 
Fig.[lk), to be later taken over by spiral waves (Fig. Q})). 
The activity time series in Fig.Hfc, displaying the fraction 
of cells where x > 0.9, shows a high frequency periodic- 
ity corresponding to the refractory period; then, the pe- 



riodicity disappears as the system settles into the spiral 
wave mode. For the time series preceding the onset of 
the spiral wave, we typically observe slow modulations 
of the envelope of the periodic oscillations, which arise 
from the interaction between the waves as well as from 
the shortcut-induced excitations. Power spectral densi- 
ties of such series were observed to show a power-law like 
decay, indicating the presence of l//-noise. 

The spiral waves were observed to be primarily created 
by a shortcut-induced excitation occurring in the refrac- 
tory "shadow" of a circular wave front, sparking a semi- 
circular wave whose transmission is partially hindered by 
the shadow. We verified that indeed spiral waves could 
be triggered in this fashion by using externally applied 
signal to stimulate an appropriate point in this region 
and then observing the resultant patterns. Once a spiral 
wave is created, it will eventually take over the dynam- 
ics of the system, as the successive excitation wavefronts 
occur with the highest frequency compared to all other 
excitations which will then be swept away. Evidently, 
the probability of a spiral wave creation per unit time 
increases with the system size N 2 and the shortcut prob- 
ability p. Thus for long times and large system sizes, exci- 
tation via random shortcuts will eventually always result 
in the formation of spiral waves. This is corroborated by 
Fig. [2j where the fraction of spiral wave configurations 
in 400 simulation runs is shown as function of time t, for 
varying system sizes at p = 0.05 (panel a) and for fixed 
system size but increasing values of p (panel b). The 
spatial structure of the spiral waves, i.e. an excitation 
front followed by a refractory shadow, makes them very 
robust against perturbations. For example, if the state 
of the system is frozen, and the state of large areas (say, 
upto a quarter of the total area) or every second cell are 
set to x = 0, the spiral wave mode is quickly recovered 
once the simulation is restarted. 

However, at high enough values of p, the shortcut- 
induced excitations become too numerous for sustaining 
the spiral wave dynamics. As almost every point is liable 
to be excited with a frequency proportional to its refrac- 
tory period, the spirals become unstable and we see a 
transition to a new regime at p = p l c ~ 0.553 (Fig. [3]). 
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FIG. 2: Fraction Fs of systems with spiral waves in 400 runs 
as function of time t: (a) fixed shortcut probability p = 0.05 
with varying system size, (b) fixed system size N 2 — 300 x 300 
with varying shortcut probability p. 
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FIG. 3: The fraction of spiral configurations Fs as a function 
of (p — pj,) normalized by system size N, for shortcut proba- 
bilities p around p l c w 0.553. Fs is calculated at t = 20, 000 
time steps, averaged over 400 runs. Inset: Dependence of the 
critical value p l c on diffusion constant D. The circles are sim- 
ulation results for N — 200 while the curve shows fitting with 
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FIG. 4: The periodic regime, recorded in a simulation of a 
iV 2 = 128 x 128 system, with shortcut probability p = 0.6. 
a)-d): Snapshots of the state of the system taken at intervals 
of At = 10 steps, e): Time series of the activity (see Fig. 1) 
and (f): the corresponding power spectrum. The main peak 
is at fo ~ 0.022, corresponding to a wavelength of « 43 time 
steps; other peaks are harmonics at integer multiples of fo. 



This value of p c was found to be independent of the sys- 
tem size TV 2 . Here, the spatial pattern becomes more 
homogeneous, as a large fraction of the system becomes 
simultaneously active and subsequently decays to a re- 
fractory state. However, some cells not participating in 
this wave of excitation carry on the activity to the next 
cycle, where it again spreads through almost the entire 
system. This results in a remarkably periodic behavior 
of the system in time, with a large fraction of cells be- 
ing recurrently active with a period close to the refrac- 
tory period of the cells (see Fig. [4] (e)-(f)). Often, small 
spiral-like waves were also observed (Fig. H (a)-(d)), but 
these were short-lived and spatial correlations could not 
be maintained. 

Finally, when p is increased still further, the very large 
number of shortcut connections guarantees almost simul- 
taneous spread of excitation to nearly all cells. As a re- 
sult, the dynamics of the system tends to "burn out" 



such that after a transient, almost all cells become re- 
fractory and not enough susceptible cells are left to sus- 
tain the excitation. We observed that once a system-size 
dependent value p — p™ (N) is approached, the probabil- 
ity of reaching the absorbing state {x 1 ^ — 0, Vi, j) grows 
rapidly (see Fig. [5|) such that at high enough values of p, 
the system is never seen to sustain its activity. The lim- 
iting behavior for N — > oo was investigated by plotting 
p™ against 1/N and extrapolating its value for 1/N — > 0. 
Fitting a quadratic function yielded p™ — > 0.86, hence it 
appears that the regime where activity is sustained has 
for all system sizes an upper limit for the shortcut den- 
sity p. Note that this can be viewed as an approximation 
only, because there is no a priori reason to assume any 
particular form for the dependence of p™ on N. 

We have also carried out simulations with disorder in 
the parameters describing the properties of individual 
cells. For instance, we have made the parameter a that 
controls the refractory period a quenched random vari- 
able ranging over a small interval. We find that there is 
a tendency for greater fragmentation of waves with dis- 
order, and corresponding increase in the number of coex- 
isting spiral waves, but otherwise no remarkable changes. 
The robustness of our results in the presence of disorder 
in the individual cellular properties underlines the rel- 
evance of this study to real-world systems, where cells 
are unlikely to have uniform properties. In addition, we 
have looked at the effect of the diffusion constant D. For 
higher D the wave propagates much faster, so that, in a 
fixed amount of time, the system will initiate excitation 
through many more long-range connections than the sys- 
tem with lower D. Therefore, the higher D system will 
be equivalent to a lower D system with a larger number 
of long-range connections, i.e., higher p (Fig. [3], inset). 

The above results carry potential relevance to all nat- 
ural systems that are excitable and have sparse, long- 
range connections. Experiments carried out in the ex- 



citable BZ reaction system with nonlocal coupling [14] 
have exhibited some of the features described above. For 
example, the coexistence of transient activity and sus- 
tained synchronized oscillations that we observe close to 
p™ has been reported in the above system, and can now 
be understood in terms of the model introduced here. 
Recently, there has also been a lot of research activity on 
the role of non-trivial network topology on brain function 
(e.g., Ref. [Hj]). I n particular, studies show a connec- 
tion bet ween the existence of small-world topology and 
epilepsy 
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161 ] . This is one of the areas where our re- 
sults can have a possible explanatory role. In the brain, 
glial cells form a matrix of regular topology with cells 
communicating between their nearest neighbors through 
calcium waves. Neurons are embedded on this regular 
structure and are capable of creating long-range links 
between spatially distant regions. It is now known that 
neurons and glial cells can communicate with each other 
through calcium waves [lTj . Therefore, the aggregate 
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FIG. 5: Inset: Fraction of configurations Fc where activity 
has ceased at t = 20, 000 time steps in 100 simulation runs 
as function of (p — Pc)- Main: The upper critical shortcut 
probability p™ as function of inverse system dimension l/N 
(o). The solid line displays a fitted quadratic function, where 
p u c (N — > oo) ~ 0.86. 

system of neurons and glial cells can be seen as a small- 
world network of excitable cells. Here the neurons, that 
are outnumbered by glial cells approximately by one or- 
der of magnitude, form the sparse, long-range connec- 
tions. Recently, the role of neural-glial communication 
in epilepsy has been investigated (e.g., Ref. [18]), but 
the observation of spiral intercellular calcium waves in 
the hippocampus [11], and the similarity of the patterns 
seen in our model with the observed features of epileptic 
seizures and bursts, makes it especially appealing to pos- 
tulate the role of small- world topology in the generation 
of epilepsy. Experimental verification of this suggested 
scenario can be performed through calcium imaging in 
glial-neuronal co-culture systems. 

The relationship between the fraction of long-range 
connections and normal brain function is indicated by 
experimental evidence not just for epilepsy but for var- 
ious other mental disorders as well. A recent study has 
found that the brains of patients with clinically diagnosed 
schizophrenia, depression or bipolar disease have lower 
glia to neuron ratio compared to normal subjects 0. 
Therefore, understanding the dynamical ramifications of 
increasing shortcuts in an excitable media can motivate 
experiments that have the potential significance of aiding 
clinical breakthroughs in treating a whole class of men- 
tal disorders. This is also connected with the question of 
the evolutionary significance for increasing glia to neuron 
ratio with brain size [20[ . It is known that excitable me- 
dia of larg er dimensions are more likely to exhibit spiral 
waves [21[. Therefore, decreasing the fraction of neurons 
(and therefore, long-range shortcuts) could be Nature's 
way of ensuring dynamical stability for neural activity. 

To conclude, we have investigated spontaneous pattern 
formation in excitable media with small-world connec- 
tions. Our results show the creation of non-trivial spatio- 
temporal patterns similar to those seen in many real-life 
systems. These patterns are formed through dynamics 



that is driven by the system's own internal architecture. 
Most important of all, the system exhibits a non-trivial 
transition point at which the pattern goes from the spa- 
tial domain, with multiple coexisting spiral waves where 
the global activity level remains more or less uniform, 
to the temporal domain, where the global activity level 
shows large oscillations as a large fraction of cells be- 
comes simultaneously active and then refractory, with a 
strict periodicity. The connection to biological phenom- 
ena, most importantly, to calcium waves in the brain, 
and the possibility of the functional role of small-world 
connections in epileptic seizures and bursts, is expected 
to make our study of special relevance. 
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